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Pseudo-CV quadratic estimators for CMB temperature and polarization power spectra have been 
used in the analysis pipelines of many CMB experiments, such as WM AP and Boomerang. In the 
polarization case, these estimators mix E and B modes, in the sense that the estimated B-mode 
power is nonzero for a noiseless CMB realization which contains only E-modes. Recently, Challinor 
& Chon showed that for moderately sized surveys {f s k y ~ 0.01), this mixing limits the gravity wave 
B-mode signal which can be detected using pseudo-CV estimators to T/S ~ 0.05. We modify the 
pseudo-CV construction, defining "pure" pseudo-C^ estimators, which do not mix E and B modes 
in this sense. We study these estimators in detail for a survey geometry similar to that which has 
been proposed for the QUIET experiment, for a variety of noise levels, and both homogeneous and 
inhomogeneous noise. For noise levels < 20 ^K-arcmin, our modification significantly improves the 
B-mode power spectrum errors obtained using pseudo-Cf estimators. In the homogeneous case, we 
compute optimal power spectrum errors using a Fisher matrix approach, and show that our pure 
pseudo-CV estimators are roughly 80% of optimal, across a wide range of noise levels. There is no 
limit, imposed by the estimators alone, to the value of T/S which can be detected. 

PACS numbers: 98.70.Vc 



I. INTRODUCTION 

In the past few years, polarization of the cosmic mi- 
crowave background anisotropy at the 10% level, a long- 
standing prediction of Big Bang cosmology 3 ha s been 
detected by several experiments 0, H3, 124 l3fj|. The 
strongest detection to date excludes zero at the 10ct level 
[32|. Going beyond first detection, a primary goal for 
CMB experiments in the next decade will be making 
precision measurements of polarization power spectra. 
In addition to providing a strong consistency check for 
the standard ACDM cosmological model, and improving 
existing uncertainties on cosmological parameters, this 
will break parameter degeneracies which cannot be re- 
solved using CMB temperature alone 01 . Examples in- 
clude constraining the reionization history of the universe 
[lii , reconstructing the dark matter distribution at 
z ~ 1 via the lensing B-mode \\M J 28 | - breaking degen- 
eracies among isocurvature modes [5| , and measuring the 
primordial fluctuations 0, ■ 

One of the most tantalizing prospects for CMB po- 
larization experiments is the possibility of detecting the 
B-mode signal from primordial gravity waves 0, l3l| . 
In the standard ACDM model, B-modes offer a unique 
observational window on these waves, since there is no 
source of gravity waves (tensor perturbations) after in- 
flation, and scalar perturbations generate only E-modes 
in linear perturbation theory. The only cosmological con- 
taminants are therefore higher-order effects such as grav- 
itational lensing of the E-mode spectrum [42^. 

On a practical level, many methods have been pro- 
posed for estimating power spectra from CMB polar- 
ization data. A fully optimal, likelihood-based anaylsis 
scales as 0(N^ ix ) and will probably not be feasible 
for next generation experiments. An alternative method, 
which has become the "industry standard" for CMB tem- 



perature experiments with N P i x > 10 5 , is using pseudo- 
Ci quadratic estimators |37|. In the polarization context, 
pseudo-C^ estimators mix E into B, in the following sense: 
for a noiseless CMB realization containing only E-modes, 
the estimated B-mode bandpowers will be nonzero. In 
an ensemble of such realizations, the estimated B-mode 
bandpowers will be zero in the mean (this is because 
pseudo-C^ estimators are unbiased by construction), but 
because they are nonzero in each realization, E-mode sig- 
nal power does contribute to the variance of the B-mode 
estimators. This extra variance can be thought of as a 
source of noise which is due to the estimators alone. Re- 
cently, Challinor & Chon {8[ showed that this "estimator 
noise" can dominate the sample variance from lensing B- 
modes, so that it becomes the dominant contaminant if 
the instrumental noise is sufficiently small. For surveys 
with f s ky ~ 0.01, they show that it limits the value of 
T/S which can be detected to ~ 0.05. 

The purpose of this paper is to construct pure pseudo- 
Ci estimators; these arc modified versions of polarization 
pseudo-C^ estimators, which do not mix E into B in the 
sense defined above. Using these estimators, the esti- 
mated B-mode power will be zero for any noiseless CMB 
realization which contains only E-modes. Strictly speak- 
ing, this is only true in the continuum limit; in a finite 
pixelization, the B-mode estimators will acquire small 
nonzero values from pixelization artifacts. We show that 
these pixelization effects can be made arbitrarily small 
by increasing the resolution. 

The basic idea of our construction was inspired by the 
pure B-mode formalism of |?ll23||. Because pure B-modes 
exist in any finite patch of sky, a BB power spectrum esti- 
mator will receive no contribution from EE signal power 
if it is built entirely out of these modes. Our pure pseudo- 
Ci estimators have the property that the observed po- 
larization map always appears contracted with a pure 
B-mode which is constructed from heuristically chosen 
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weight functions. The main technical difficulty is en- 
suring that the fast algorithm for calculating the trans- 
fer matrix [T2I Appendix E] still goes through after this 
modification. 

Throughout this paper, we use a fiducial ACDM cos- 
mology which is consistent with WMAP [33|, with pa- 
rameters { fifc/i 2 , f2 m /i 2 f^A, t, |Ai?| 2 , n s , w, vn v } = { 
0.024, 0.14, 0.73, 0.17, 2.57 x 10~ 9 , 1, -1, }. For pix- 
elized all-sky maps, we use the Healpix coordinate system 
[Til HI ] exclusively. We use a phenomonological defini- 
tion of T/S, defining it to be the ratio of temperature 
multipoles C|f = n ^ or /C|^ or . 

This paper is organized as follows. In WI\ we briefly 
review pseudo-CV estimators, before constructing pure 
pseudo-Cf estimators in 9IVI In we consider a spher- 
ical cap shaped mock survey with uniform white noise, 
using it to illustrate general features of our estimators, 
and show that our modification significantly improves the 
performance of the estimators for noise levels < 20 /iK- 
arcmin. As a first step toward more realistic instrumen- 
tal noise, in WII we consider a survey with inhomoge- 
neous, but not spatially correlated, noise and show that 
the same conclusions apply. In 3VIII we study a range 
of noise levels and show that our estimators are roughly 
80% of optimal, defined by the degradation in the value 
of (T/S) which can be detected, for all noise levels con- 
sidered. We conclude in Willi 



E and B modes are defined as follows. First, we define 
operators 

1 



Sab = -VaVft + ^ffafcV 2 
B ab = ^£acV C V fc + ie fcc V C V a 
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which take scalar fields to symmetric traceless tensors. 
E and B harmonics are defined by 
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We will also need some analagous definitions for the spin- 
1 case. The "two-point" basis X a , Y a has already been 
defined; a global basis x a , y a is defined by: 
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y a = (0 sin0). (5) 



Spin-1 harmonics, which we label G and C for "gradient" 
and "curl" , are defined by 
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II. NOTATION AND CONVENTIONS 

We represent CMB polarization as a symmetric trace- 
less tensor IT a (, as in Kamionkowski, Kosowsky fc Steb- 
bins |19| , but we have changed some conventions to agree 
with those of CMBFAST and Healpix. The all-sky metric 
g ab and antisymmetric tensor e ab are given by 
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We define basis polarization fields 
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Given two points x, x', we define X a to be the vector at x 
which points away from x' along the great circle arc which 
connects the two, vector Y a = —e a f,X b , and symmetric 
traceless tensors Q a b = (X a X b -Y a Yb)/2, U a b = (X a Y b + 
Y a X b )/2. The quantities X', Y', Q', U' are defined in the 
same way with x, x' interchanged. Note that we always 
use lower case to distinguish the basis q ao , u ab , which is 
globally defined (except at poles), from the "two-point" 
basis Q a b, U ab , which is defined relative to a pair of points 
x. x'. 



III. PSEUDO-Cf ESTIMATORS 

Pseudo-Cf quadratic estimators for polarization power 
spectra have been previously constructed by Hansen & 
Gorski [l~2j], who used the pseudo-CV formalism In 
this section, we briefly recall the pseudo-CV construction, 
largely for the sake of establishing notation. 

One first chooses a pixel-dependent weight function 
W(x) which is zero outside the survey region. The choice 
is made heuristically, but the performance of the estima- 
tors is improved (on angular scales with signal-to-noise 
ratio < 1) by choosing VT(x) to be smaller where the 
noise is larger, in order to downweight noisy regions. Fre- 
quently, the weight function is also apodized near the sur- 
vey boundary, in order to reduce harmonic ringing. One 
then defines pseudo E and B multipoles by (the tildes 
signify "pseudo") 

E lm = 5>n ab (x)iF(x)y ( £,; )ab (x) (8) 

X 

B lm = ^2H ab (x)IF(x)y ( f r : )Qb (x). 

X 

and pseudo power spectra by 
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It can be shown that the expectation values of the pseudo 
power spectra are given by 
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are 
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max-by-f-max transfer matrices and Nf E , 
-by-1 vectors which represent additive 
noise bias. There is an efficient algorithm, which will 
be discussed in Appendix [D] for exactly computing the 
transfer matrices from the weight function. The noise 
bias can be computed exactly in cases where the noise is 
uncorrelated between pixels (Appendix [EJ) , or by Monte 
Carlo in the general case. 



Unbiased power spectrum estimators C 1 
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be obtained from the pseudo power spectra Cf E , Cf B 
by simply subtracting the noise bias and applying the 
inverse of the (2£i max )-by-(2£ max ) "grand unified transfer 
matrix" : 
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The preceding construction has assumed that the 
power spectrum is estimated at every multipole I. For 
reasons of sky coverage or signal-to-noise, it is often nec- 
essary to bin multipoles into bandpowers with A£ > 1. In 
this case, for each band b, one defines pseudo bandpowers 



p^tEE \ r, p^EE 



CfB^PUC?*, (12) 



where the matrix P defines the £ weighting within each 
bandpower estimator. A commonly-used choice is [15J: 
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One also introduces a matrix P, which defines an interpo- 
lation scheme by which the signal power spectra depend 
on bandpowers A/,: 
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and are related to the unbinned versions by: K bb , = 

PbiK u ,P£'b', Nb = PuNi. The binned analogs of the 
unbiased estimators (Eq. ((TTJl ) arc defined by: 
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Strictly speaking, these are unbiased estimators of the 
bandpowers Ah, when the CMB power spectra are of the 
precise form Q14[l. 



IV. CONSTRUCTING PURE PSEUDO-C^ 
ESTIMATORS 

We now construct a modified version of the B-mode 
power spectrum estimator which is "pure" , in the sense 
that the estimator is identically zero in any noiseless re- 
alization of the CMB which contains only E modes. It 
is convenient to explain the construction first in the con- 
tinuum limit ( 8IVA| . where sums over pixels can be re- 
placed by integrals, and then address implementational 
issues associated with finite pixelized maps f Ul V B|l . 



A. Constructing pure pseudo-Cj estimators 1: 
continuum limit 

The underlying reason why the unmodified pseudo- Ci 
estimator is not "pure", in the sense defined above, can 
be understood as follows. The pseudo B multipoles are 
defined by 



B tm = J d 2 x^2n ab (x)W(x)Y B ^ )ab (a 



(18) 



Since multiplication in position space mixes E and B 
modes, the product W(x)Y B m - )ab (x) is a mixture of E 

and B, even though Yam)ab^ x ) * s a B-mode. Therefore, 

Bg m receives a nonzero contribution from E-mode power 

in n ah . 

The basic idea of our construction is now easy to ex- 
plain. Writing out the definition of Y B m ^ from Eq. 

the definition of Bi m above can be rewritten 



A commonly-used choice, corresponding to piecewise flat 
power spectra, is: 
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lb 
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1(1+1) 11 how - * ^ How (15) 
otherwise. 



The binned analogs K bb , , Nh of the transfer matrices and 
noise bias vectors (Eq. Q1U[0 arc defined by: 
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Suppose that this definition is modified by simply moving 
W(x) inside the operator B a b- 



j^pure def 
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d 2 x ^2U ab (x)B ab (W(x)Y; m (x)) 



(20) 
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In addition, suppose that the weight function W{x) and 
its gradient vanish on the boundary of the region. Then 
B^ e will be identically zero in a noiseless CMB real- 
ization which contains only E modes. One way to see 
this is from the perspective of the pure B-mode formal- 
ism 0,1131; since (W(x)Yi m {x)) satisfies both Dirichlet 
and Neumann boundary conditions, B a b{W(x)Yi m {x)) is 
a pure B-mode. Alternately, a short, self-contained proof 
is given as follows. A noiseless CMB realization which 
contains only E-modes can be written W a b — B ao <\>, for 
some scalar function <f)(x). Putting this into the defini- 
tion of and integrating by parts twice, one gets 



^(l- 1)1(1 +!)(* + 2) ~ pure 



J d 2 Xy/g 
J d 2 xyfg 
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= 0. 



Here, §da denotes integration around the survey bound- 
ary and t a denotes the unit tangent vector. The bound- 
ary terms are zero because W{x) and V a W(x) are as- 
sumed to vanish on the boundary. 

To avoid confusion, we will denote the "pure" pseudo 
B multipole defined in Eq. (|2*U1) by Bf^[ e , and denote 
the "mixed" multipole defined by Eq. (JTSJ by B%£ xed , 
for the rest of this paper. This one change, replacing 
j^mixed B^ e , is the only modification we propose to 



the pseudo-C^ formalism of £11111 In particular, we leave 
the definitions of Eg m and Cf E unchanged. 

By analogy with the definition of Bf^ e , it is possi- 
ble to define a "pure E" multipole E^ e , but we expect 
that incorporating this would worsen the performance of 
the EE power spectrum estimator. This is because am- 
biguous modes, which are useful for measuring E-mode 
power, would be discarded. Note that we have named 
our estimators "pure pseudo-CV" even though we are re- 
ally using a "mixed" estimator for E-modes and a "pure" 
estimator for B-modes. 



B. Constructing pure pseudo-CV estimators 2: 
finite pixelized maps 

Before moving on to construct pseudo power spectra 
and transfer matrices for the modified estimators, we ex- 
plain how the definition of B^ e (Eq. Q2U[l) is to be in- 
terpreted when Tl ab (x) is a finite pixelized map, instead 
of an idealized continuous field. The technical obstacle 
is making sense out of the object B a b(WYg m ) which con- 
tains covariant derivatives, and appears on the right hand 
side of Eq. J30J. We first expand B a b(WYi m ) using the 
product rule 

B ab {hh) = {Babh)f2 + fi{B ab f 2 ) (22) 
+2T a6cd (V7i)(V d / 2 ), 

where we have introduced the tensor 

Tabcd — -^acQbd + ~^ad9bc + ^bc9ad + ~^bd9ac- (23) 

(Note that contracting a spin-2 object with T results in 
a 45° rotation, e.g. T abc dq cd = -u a b, T abc d£ cd = ~B a b-) 
Performing the derivatives which act on Y^ m , and re- 
placing the integral by a sum over pixels, this results 
in the following form for B^ e . (We do not include the 
pixel area as an overall prefactor, since the normalization 
of B^" e will eventually drop out when we define unbi- 
ased power spectrum estimators in Eq. H3U|) . However, 
straightforward conversion of the integral to a sum does 
implicitly assume an equal- area pixelization.) 



(24) 
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z ^2T»^H o6 (x)W c (x)y ( ^ )d (x) 
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where we have defined 
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W a = V a W 

W ab = (V Q V 6 -(l/2) 5ah V 2 )t¥ 

N' t = 2/ V(l-l)(l + 2) 

N e = l/v/(<-W+l)(< + 2) 



In the form l|24|) , the only covariant derivatives are those 
which act on the weight function in the definitions Ij25(l 
of W a , W a b- If W(x) is of known analytical form (e.g., 
Gaussian or cosine apodization), then W a , W ab can sim- 
ply be computed analytically. More generally, assuming 
that W(x) is slowly varying compared to the pixel scale 
(which will be the case for any sensible weighting), W a 
and Wab can be computed from W(x) by finite differ- 
encing. This is nontrival to implement for an irregular 
spherical pixelization such as Healpix; we present one fi- 
nite differencing scheme in detail in Appendix iBl 
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In either case, once a prescription for computing W a , 
W a t from W has been specified, we use Eq. 124(1 as the 
definition of B v e ^ e in a finite pixelization. Note that the 
first term on the right-hand side is the unmodified pseudo 
multipole B^ J l l xed ; the second and third terms can be 
thought of as counterterms which cancel the E-B mixing 
and involve spin-1 and spin-2 weight functions W a , W a b- 
The relative weights of the three terms are ^"-dependent, 
in a way which downweights the counterterms at high I. 
This is in qualitative agreement with the general picture 
of E-B mixing arising from survey boundaries 6, 23] in 
which the mixing is a less significant contaminant at high 
L 

The definition (|24() suggests an efficient algorithm for 
evaluating the pseudo multipoles -Bf" re . The first term 
on the right-hand side is computed, for all (£, m) simulta- 
neously, using a fast spin-2 spherical harmonic transform 
of the weighted polarization field W(x)H a b(x). The sec- 
ond and third terms are analagously computed using a 
spin-1 transform of the vector field 



pseudo spectrum Cf E are given by 



2 r T a Yl a0 W c = (UqW y - UuW x )x d 

+ (Tl Q Wx+Tl u W Y )y a 



(26) 



and a spin-0 (ordinary spherical harmonic) transform of 
the scalar function 



2T abcd X\ ab W cd = TIqWu - TIuWq. 



(27) 



(In Eqs. (|26|l and (|27l) . we use subscripts X, Y to denote 
components of a vector field in the "global" {x a ,y a } ba- 
sis, and subscripts Q, U to denote components of a spin-2 
field in the {q a b,u a b} basis.) The prefactors N' e , Ni are 
applied after the transforms. This algorithm computes 
every Bf^ e from the map Tl a i,(x), with the same asymp- 
totic complexity C(^m aa; ) as in the unmodified pseudo-CV 
formalism of £11111 but the overall constant is 2-3 times 
worse, since one spin-0, one spin-1, and one spin-2 har- 
monic transform are needed, rather than a single spin- 
2 transform. Since fast spin-1 harmonic transforms are 
not implemented in the Healpix library, we implemented 
them as part of this paper; the details are presented in 
Appendix [X] 

We have now constructed "pure" pseudo multipoles 
Bf^ e for finite pixelized maps; we conclude this section 
by defining pure pseudo power spectra and unbiased esti- 
mators of the power spectrum. This is done in complete 
analogy with flllll The pseudo power spectrum (j BB 'P ure 
is defined by 
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The expectation values of (j BB 'P ure anc j ^ ne unmodified 
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An essential ingredient in the method is a fast algorithm 
for computing the new transfer matrices K^ g f ure . In Ap- 
pendix^ we present such an algorithm and show that 
the computational cost is 0(^ aa .), which is the same as 
the cost of evaluating the estimators once. (Noise bias 
is discussed in Appendix IE!) The matrix K^ e f ure mea- 
sures the contribution to {Cf B,pure ) by E-mode power, 
and is therefore expected to be zero. Strictly speaking, 
this is only true in the continuum limit; in a finite pix- 
elization, K^ e f ure will acquire a small nonzero value from 
pixelization artifacts. This will be studied quantitatively 
in El 

Finally, we define pure pseudo-CV estimators (j EE 'P ure ) 
(jBB,pure ^ s i m ply subtracting noise bias and applying 
the inverse of the (2£ max )-by-(2£ max ) transfer matrix: 



EE. pure 
i 

p^BB,pure 



K + K7, 

r^—pure r^-\-pure 



(30) 



Cf, B - 



-N EE 

jyBB ,pure 



Note that, even though we have not changed the def- 
inition of the pseudo spectrum Cf E , the definition of 
the unbiased estimator Cf E has been modified in our 
formalism. This is because the matrix inversion on the 
right-hand side of Eq. (|30jl mixes all rows of the (2£ max )- 
by-(2£ max ) matrix. In practice, we have found that the 
change in Cf E is miniscule, but we do introduce the no- 
tation (j EE -P ure to distinguish between the two versions. 



V. AN EXAMPLE WITH HOMOGENEOUS 
NOISE 

In the preceding section, pure pseudo-CV estimators 

£,EE,pure^ gBB.pure ^ w JjJ cn nQ ^ m j x £ &n( \ g modes 

have been defined. In this section, these estimators will 
be studied in detail for a specific mock survey, namely 
a spherical "cap" of radius r = 13°, with uniform white 
noise, and a Gaussian beam with 25' FWHM. Parameters 
similar to these were originally proposed for the QUIET 
experiment [45]. which provided the original motivation 
for this paper. (We note that since completion of this 
paper, the proposed QUIET beam size has changed to 
10'.) Power spectra will be estimated in "flat" (Ci cx 
!/£(£ +1)) bands with Al — 40. An exception is the 
lowest band, where we have restricted the £ range to 10 < 
£ < 40; we have found that trying to estimate power at 
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I < 10, which is above the survey scale, gives spurious 
results. 

We will assess the performance of the estimators by 
comparing their Monte Carlo covariance with the Fisher 
matrix, 

F w = ^(SfcCSo + N^Sb^So+N)- 1 ) (31) 

Here, N is the (2N p i X )-by-(2N p i X ) noise covariance ma- 
trix, So is the signal covariance matrix in the fiducial 
model, and S;, is the signal covariance matrix associated 
to each fiat bandpower in EE or BB. The Cramer-Rao 
inequality asserts that the Fisher matrix is a lower bound 
on the covariance of any unbiased power spectrum esti- 
mator; conversely, minimum- variance unbiased quadratic 
estimators [34| have a covariance matrix which is equal 
to the Fisher matrix, but require dense matrix compu- 
tations which are prohibitively expensive for the survey 
sizes considered in this paper. 

The infeasibility of dense matrix algebra also means 
that the Fisher matrix (|3 lf> cannot be computed in a 
straightforward fashion. Our method for making the cal- 
culation affordable is to exploit azimuthal symmetry of 
the survey region and noise. Taking a Fourier transform 
in the azimuthal coordinate tp, the matrices S, N will 
be block diagonal in the azimuthal wavenumber m (S is 
still dense in the coordinate 9), which makes the matrix 
operations in Eq. (|31[1 affordable. Details of the method 
are presented in Appendix [FJ We note that this method 
also could be used to make optimal power spectrum es- 
timation affordable, but it applies only for surveys in 
which both the sky coverage and noise are azimuthally 
symmetric. Such surveys therefore permit benchmark 
comparisons between pseudo-CV estimators and optimal, 
likelihood-based methods, for survey sizes which are large 
enough that comparison would normally be infeasible. 

To use the pseudo-C^ method, one must heuristically 
choose a weight function W(x). We will make a modest 
effort to optimize this by hand, deferring a systematic 
framework for a future paper. The only possibility con- 
sidered for W(x) will be cosine apodization, 



1 



for 9 < r — 



W(0,<p) 



i - i cos U^J for r - r* < 6 < r 
for 9 > r 



The parameter r* is an apodization length which will be 
chosen shortly. Note that this choice of W(x) satisfies 
the requirement of EHVI both W and its gradient vanish 
at the survey boundary. We obtain the spin-1 and spin-2 
weights W a and W a b (Eq. by analytic calculation, 

rather than finite differencing: 



W a = - — am | 7T- 



W ab = 



7r r 
23 - S |, 



7rcot(0) / r 
' ■ sm 71" — 



2r* 



(33) 



for r — r* < 9 < r, and otherwise. 

We note that the "pure" multipole B^ e is defined 
(Eq. EU as the overlap integral between U ab (x) and the 
pure B-mode 

W(x)Y { f m)ab (x) + N' e T ab cd W c (x)Y^ m)d (x) (34) 
+N e T ab cd W cd (x)Y lm (x). 

It is illuminating to consider the behavior of this mode 
for varying £ and r* . In Figure ^ the mode is shown for 
(£,m) = (10,1) and (£,m) = (50,1), taking r* = 7°. At 
i = 10, the second and third terms in l|34() dominate the 
first, and the mode is concentrated in the apodization 
region (r — r*) < 9 < r, which is undesirable from a 
signal-to-noise perspective. At £ = 50, the contribution 
of these terms is comparable to the main term, and the 
statistical weight is distributed throughout the survey 
region. 

This illustrates a general point: the "pure" estima- 
tors defined in fllVI require apodization near the survey 
boundary. In the tophat limit — > 0, the derivatives of 
W(x) increase in magnitude; the third term in the mode 
(|34|l dominates the others, and is concentrated at the 
boundary. In this limit, formally receives no con- 

tributions from E-modes, but achieves this using coun- 
terterms which are line integrals around the boundary, to 
cancel the E-B mixing. (Compare Eq. (4) of Lewis |2l|-1 
Naturally, for noisy data, this ruins the performance of 
the power spectrum estimator. In a finite apodization, 
the counterterms in l|34|l are effectively "smeared" over a 
nonzero area near the boundary, and the resulting power 
spectrum estimators are sensible in the presence of noise. 
By comparison, unmodified pseudo-CV estimators behave 
reasonably even for a tophat weight function; apodization 
serves the milder purpose of reducing Fourier ringing. 

Another general feature of our method is that more 
apodization is optimal at low £. At £ = 10, Figure 
shows that the counterterms in l|34|) are dominant, if r» = 
5°. Since increasing decreases the magnitude of the 
counterterms, this suggests that the optimal apodization 
length at I = 10 is greater than r» =5°. In contrast, at 
£ = 50, FigureQ]suggests that r» ~ 5° is roughly optimal. 

Next we show the behavior of the transfer matrices 
K^,, K bb f ure ■ In Figure [3 the matrix entries are shown 
for varying with the band b fixed at £ m i n = 80, 
£max = 120. These matrix entries can be interpreted as 
the contribution to (Cjf E ) and (fj^ B 'P ure ^ f rom E-mode 
and B-mode power in the band 6'. The small nonzero 
value of K b y Ure represents contamination of the B-mode 
estimators by E-mode power in a finite pixelization, and 
goes to zero in the continuum limit N S id e — ► oo. 

The transfer matrix formalism can also be used to cal- 
culate the total contribution to ((j BB 'P ure ^ f rom & \\ e_ 
mode power in the fiducial model, including power at 
£ > 200. One can think of this as the pseudo power 
spectrum of E-mode power which is aliased to B by the 
pixelization. The result is shown in Figure|3 for compar- 
ison, a gravity wave B-mode spectrum with T/S = 0.01 



FIG. 1: The pure B-mode 13411 . for (£,m) = (10, 1) (left) and (50, 1) (right), using apodization length r* =5°. Only the real 
part is shown. In the left panel, the second and third terms, which are concentrated in the apodization region, dominate the 
main term. 



1 

0.1 
0.01 
0.001 
0.0001 

10- 5 b 

10- 6 r 
IO" 7 - 



K 



bb' 







K 



bb" 



50 



1 

0.1 
0.01 
0.001 
0.0001 

10- 5 
io- 6 
io- 7 



K 



pure 
bb' 



K b -r e (N slde =i2s; 



(N slde =512) - 

J I I I I I I I I I I I I I L 



100 

Multipole 1 



150 



200 50 100 150 200 

Multipole 1 



FIG. 2: Transfer functions Kj~i, (left) and K,y UrR (right), shown for varying b', with the band b fixed at l m in = 80, £ m ax = 120. 
The apodization length used was r* = 5°. The matrix entry K^ b , represents the mean response of the BB estimator in band b 
to BB power in band b'\ Kr,, represents the response of the BB estimator to EE power. As N si d c —* oo, KZ^ ure approaches 
zero. 



is also shown. As expected, the aliased power goes to zero 
in the continuum limit. The level of aliased power, rela- 
tive to the T/S = 0.01 gravity wave signal, suggests that 
detecting such a small B-mode signal requires a massive 
amount of overpixelization. We will quantify this better 
in Will 

We now take up the issue of choosing the apodization 
length r* in Eq. Ip^, In Figure H the Monte Carlo RMS 
scatter of the unbiased estimator (j BB ^ ure j s shown for 
varying r», in two bands b: the lowest band 10 < £ < 40, 
and the second lowest band 40 < £ < 80. For both bands, 
the estimator performance degrades sharply when is 



chosen smaller than the optimal value, but the optimal 
value is r* ~ 8° for the first band and ~ 3° for the 
second. This is consistent with the qualitative discussion 
after Eq. (|34|l . For the higher bands, repeating this anal- 
ysis shows that r* ~ 3° is roughly optimal. Our solution 
to this problem is to use one weight function W(x), with 
r* = 8° for the lowest i band, and a different W{x), with 
r* = 3° for the higher bands. 

We briefly describe the implementational details asso- 
ciated with using different weight functions in different 
I bands. We emphasize that these comments apply to 
all types of pseudo-C^ estimators, not just the pure es- 
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FIG. 3: Contribution to (C b ,jmre ) from the fiducial E-mode power spectrum, for zero beam (left) and FWHM 25' (right). 
This can be interpreted as the estimated B-mode power which is contributed by E-modes due to pixelization artifacts. A 
gravity wave B-mode spectrum with T/S = 0.01 is shown for scale. The apodization length used was r* = 5°. 
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FIG. 4: RMS bandpower errors for pure pseudo-CV estimators (blue/solid), unmodified psuedo-Cf estimators (green/dashed), 
and optimal Fisher errors (red/dotted), shown for varying apodization length r*. In the left panel, the band 10 < £ < 40 is 
shown; in this band, the performance of unmodified pseudo-C^ estimators is very poor and the dashed curve is not visible. In 
the right panel, the band 40 < £ < 80 is shown. The noise level used was 20 pK-arcmin. 
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timators studied in this paper. To evaluate the estima- 
tors, one first computes a separate set of Bi m for each 
weight function. Then, in each band b, the pseudo power 
spectrum Cb is computed using the Bi m corresponding to 
the appropriate weight function. Generally speaking, one 
could use any number N w t < Nband of weight functions, 
at the expense of increasing the computational cost, since 
N wt sets of spherical harmonic transforms must be com- 
puted. (The increase in the cost is smaller than a factor 
N wt , since the transforms are 0(£? max ). For the specific 
example in this section, using a different weight function 
in the lowest band, with £ max — 40, has a very small im- 
pact on the running time.) When computing the transfer 
matrix element Kbb', one uses the weight function asso- 
ciated to band b. 

We have now specified the weight functions which 
will be used for the mock survey in this section. To 
make a fair comparison between the pure and unmodi- 
fied pseudo-Cf estimators, we optimized the apodization 
length r* independently for the unmodified versions, but 
found that a tophat weight function (r* = 0) was opti- 
mal in all bands. In the lowest £ band (10 < £ < 40), we 
found that unmodified pseudo-CV estimator performance 
was extremely poor (Figure^) . If E- modes are artificially 
removed from the signal, then it improves dramatically, 
suggesting that this is due to severe E-B mixing in the 
lowest band. 

In Figure [SJ we have compared (Monte Carlo) RMS 
errors of pure pseudo-Q estimators, RMS errors of the 
unmodified pseudo-C^ estimators, and the Fisher bound 
(|31[) . This was done for two choices of signal and noise 
level: first, a T/S = 0.2 gravity wave signal and noise 
level 20 /xK-arcmin, and second, a T/S — 0.05 gravity 
wave signal and noise level 10 /xK-arcmin. For E-modes, 
the performance of the pure and unmodified estimators 
is the same, as stated at the end of flIVI and close to 
optimal except in the lowest band. For B-modes, the 
relative performance depends on the noise level: for the 
first choice, the pure estimators perform slightly better 
than the unmodified versions, and slightly worse than 
the Fisher bound, in all bands. The exception is the 
lowest band, in which the performance of the unmodified 
estimator is very poor. For the second choice of noise 
level, the performance of the unmodified estimators has 
degraded significantly, but the performance of the pure 
estimators remains about the same relative to optimal. 



VI. AN EXAMPLE WITH INHOMOGENEOUS 
NOISE 

In the previous section, we studied a mock survey with 
homogeneous white noise. As a first step toward more re- 
alistic noise models, in this section we study noise which 
is inhomogeneous, but not correlated between pixels. We 
will use the same region as in Jvl namely a spherical cap 
of radius r = 13°, but the following special form for the 



noise covariance: 



(Q(x)Q(x')) = (C/(x)C/(x')) 



n 



(35) 



1 — cos r 



sin(0) <5 2 (x-x') 



Here, 77 is a constant with units /iK-arcmin which we 
will use to quote the noise level. The normalization 
r/(l — cosr) is included so that the total integrated sen- 
sitivity will be the same as a homogeneous survey whose 
noise level is 7/. As in fJ3 we will study two sensitiv- 
ity levels: 7/ = 20 /iK-arcmin and r\ = 10 /xK-arcmin. 
Because the noise covariance (|35l) is azimuthally sym- 
metric, optimal power spectrum errors can be calculated 
using the method of Appendix iFl and used to benchmark 
pseudo-C^ power spectrum estimation. 

To motivate the form (|35|l . consider a single detector 
which makes constant-velocity scans through the center 
of the survey region, at a variety of angles. If the distri- 
bution of angles is uniform, the scans are rapid, and the 
timestream noise is white, then the resulting noise in the 
pixel domain will be given by l)35|l . 

We must first choose pixel weight functions which will 
be used to construct pseudo-CV estimators for the noise 
(|35p. Considering E-mode estimators first, we have found 
that a tophat weighting is significantly suboptimal at 
high £, in contrast to the homogeneous case. (A sim- 
ilar result for temperature estimators appears in [l3j . 
where it is shown that Gaussian apodization improves a 
tophat weight function in a mock survey with inhomoge- 
neous noise, but not in the homogeneous case.) However, 
we were able to construct E-mode estimators which were 
95% of optimal, for all £ > 80, by using weight functions 
of the form 



W e (6) 



sin(0) 



(36) 



The form of W 6 {&) was motivated by the observation that 
inverse noise weighting, W(9) — 1/ sin(0), should be opti- 
mal in the noise-dominated limit; we included a regulator 
e to smooth the singularity at 9 — 0. The values of e used 
are shown in Tabled We were unable to achieve the same 
level of EE estimator performance using cosine apodiza- 
tion (|32|l in place of the apodization (|36[) . or using fewer 
than four different weight functions. 

This example highlights a practical issue for pseudo-CV 
power spectrum estimation: for a given survey, how does 
one choose pixel weight functions which optimize perfor- 
mance of the estimators? This issue is separate from the 
E-B separation problem which is the focus of this paper; 
we have seen here that it is important even for E-mode 
power spectrum estimation alone. In a future paper, we 
plan to investigate algorithms for constructing optimal 
pixel weight functions, starting from Monte Carlo simu- 
lations of the noise. 

Turning now to weight functions for B-mode estima- 
tors, we have found that cosine apodization i|32|) re- 
sults in pseudo-Cf estimators whose performance (rela- 
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FIG. 5: RMS bandpower errors for pure pseudo-CV estimators (blue/left), unmodified pseudo-CV estimators (green/middle), 
and optimal Fisher errors (red/right). The survey region is a spherical cap of radius 13° with uniform white noise. The top 
two panels are for a fiducial model with T/S = 0.2 and noise level 20 /iK-arcmin; the bottom two panels are for a fiducial 
model with T/S = 0.05 and noise level 10 /iK-arcmin. In each pair of panels, bandpowers are shown for E-modes (left) and 
B-modes (right). In the E-mode panels, multipoles £ < 250 have been replotted on a log scale for visibility. In the B-mode 
panels, the lensing component of the power spectrum has been shown separately (dashed); detecting the gravity wave signal 
requires measuring power in excess of this level. 



tive to optimal) is similar to the homogeneous case. How- 
ever, best results were obtained using apodization length 
r* = r = 13° in all bands, in contrast to the homoge- 
neous case where we used r* = 8° for 10 < £ < 40 and 
r* = 3° for i > 40. 



We have now completely specified weight functions; 
the performance of the pseudo-CV estimators is shown 
in Figure The results are similar to the homoge- 
neous case. At sensitivity r\ — 20 /iK-arcmin, the pure 
B-mode estimators perform slightly better than the un- 
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rj = 20 /iK-arcmin 
£ range weight function 



£ < 440 tophat weighting 
400 < £ < 560 W e {6), e = 0.070 
560 < I < 680 W e (6), e = 0.025 

^>680 W e (9), e = 0.010 



r\ = 10 /iK-arcmin 



£ < 600 tophat weighting 
600 < £ < 720 W e {6), e = 0.070 
720 < £ < 800 W e (6), e = 0.017 

l> 800 We(0), e = 0.008 



TABLE I: Pixel weight functions used to construct E-mode 
estimators for the mock surveys in this section. 



modified estimators, and slightly worse then optimal. At 
increased sensitivity r\ = 10 /iK-arcmin, the performance 
of the pure B-mode estimators remains the same relative 
to optimal, but the unmodified estimators have signifi- 
cantly degraded. Our conclusion in this section is that 
inhomogeneous noise, at least for the specific form 1|35|) 
studied here, does not pose a problem in principle for 
power spectrum estimation using pure pseudo-CV esti- 
mators, although the practical issue of choosing suitable 
pixel weight functions remains. 



VII. DETECTABILITY OF GRAVITY WAVE 
B-MODES 

In this section, we will study the smallest gravity wave 
B-mode signal which can be detected using our pure 
pseudo-Cf estimators, as a function of the noise level. 
We use the survey parameters of ^3 a 13° spherical cap 
with uniform noise, and define pseudo-C^ estimators us- 
ing the I bands and pixel weight functions defined in that 
section. 

In the left panel of Figure [7J we have considered a 
"wide-beam" survey with a 25' Gaussian beam. The 
lcnsing B-mode signal is treated as an extra source of 
Gaussian noise. The minimum (T/ S) which is detectable 
at la is given by 



(T/S) la = (5E) b C^(5Ef 



-1/2 



(37) 



where C bb denotes the bandpower covariance of the esti- 
mators (including noise contribution) in a fiducial model 
with T/S = 0, and (5E) b denotes the estimator mean 
contributed by a tensor S-mode signal with T/S = 1. 
We note that (SE) b can be computed exactly using the 
transfer matrix formalism of Appendix [Dj but C bb must 
be computed by Monte Carlo. 

In Figure \7\ (left panel) , we show the minimum de- 
tectable T/S for both pure and unmodified pseudo-CV 



estimators, with the optimal Fisher value (computed us- 
ing the method of Appendix [FJ) shown for comparison. 
In all three cases, there is a floor to the gravity wave sig- 
nal which can be detected, arising from lensing contami- 
nation, even in the limit of zero noise. With unmodified 
pseudo-C^ estimators, we find the best possible la detec- 
tion is T/S = 0.042, which is reached at a total sensitivity 
of ~ 10 /iK-arcmin. This agrees well with the results of 
Challinor & Chon With pure psuedo-C^ estimators, 
the best possible la detection is T/S = 0.00170, which is 
reached at total sensitivity ~ 2 /iK-arcmin. For compar- 
ison, the optimal value is T/S = 0.00133, or 78%, so the 
performance of the estimators is close to optimal. These 
(T/S) values are specific to the survey region considered 
here, and depend on the geometry in a way which in- 

— 1/2 

volves boundary E-B mixing and does not scale as f sky 
as mode-counting arguments would suggest 0. 

In principle, these values of (T/S) are not necessarily 
ultimate limits, since we have arrived at them by treating 
the lensing B-modes as a Gaussian contaminant which 
can not be separated from the gravity wave signal. In 
actuality, lensing B-modes are non-Gaussian, and "de- 
lensing" algorithms have been proposed 0, [28| which 
exploit this to separate the gravity wave and lensing 
components, reducing the level of contamination. The 
amount of lensing contamination which can be removed 
in this way depends on the noise level and beam size. In 
this paper, we make no attempt to model residual dclcns- 
ing errors realistically. Instead, we note that the result 
of delensing will lie between two extremes: no separation 
of the gravity wave and lensing B-modes (which has just 
been considered), and complete separation. In the right 
panel of Figure we consider the latter extreme; we 
completely remove the lensing component of the B-mode 
power spectrum, and show (T / S)\ a as a function of noise 
level. We have also used zero beam size, to crudely re- 
flect the fact that delensing algorithms require measuring 
modes at high £ to separate the gravity wave and lensing 
signals at low I. Under these extreme assumptions, there 
is no limit in our framework to the detectability of grav- 
ity wave B-modes as the signal-to-noise improves, since 
we do not consider other sources of contamination such 
as astrophysical foregrounds. 

Using unmodified pseudo-C^ estimators, the smallest 
gravity wave signal which can be detected is the same 
as in the left panel, T/S = 0.042; removing the lensing 
component of the B-modes results in no improvement. 
Rephrasing, the extra BB estimator covariance which 
is contributed by E-modes dominates the contribution 
from lensing B-modes. Using pure pseudo-C^ estimators, 
there is a floor to (T / S)\ a in any fixed pixelization, but 
the floor improves without limit as the resolution is in- 
creased (N S ide — 256,512, 1024 are shown.) This is con- 
sistent with the discussion in Jvl in a finite pixelization 
there is some contamination of the B-mode estimators by 
E-modes, but the contamination goes to zero in the con- 
tinuum limit. For any fixed noise level, one can choose a 
sufficiently high resolution so that the contamination is 
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FIG. 6: RMS bandpower errors for pure pseudo-CV estimators (blue/left), unmodified pseudo-Cf estimators (green/middle), 
and optimal Fisher errors (red/right). The survey region is a spherical cap of radius 13° with inhomogeneous noise given by 
1351 1. The top two panels are for a fiducial model with T/S = 0.2 and noise level 20 /^K-arcmin; the bottom two panels are for 
a fiducial model with T/S = 0.05 and noise level 10 /^K-arcmin. In each pair of panels, bandpowers are shown for E-modes 
(left) and B- modes (right). 



negligible. If this is done, then the value of (T/S)i a ob- 
tained is approximately 80% of the optimal Fisher value, 
for all noise levels considered in the right panel of Figure 

13 

The results of this section illustrate a qualitative dif- 
ference between our pure pseudo-C^ estimators and the 
unmodified versions. As the signal-to-noise improves in a 



fixed survey region, a floor is revealed to the gravity wave 
signal which can be detected using unmodified pseudo-C^ 
estimators. This is not the case for pure pseudo-CVs, al- 
though a second smaller floor is eventually revealed when 
the signal-to-noise becomes good enough that the sensi- 
tivity is limited by contamination from lcnsing. Going 
beyond this will require separating the lcnsing signal us- 
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FIG. 7: Minimum T/S detectable at la, for a 13° spherical cap survey with uniform noise, assuming a 25' beam and Gaussian 
lensing contaminant (left), and zero beam and no lensing contamination (right). The solid line is the optimal Fisher value, 
the dotted line is the Monte Carlo value using unmodified pseudo-CV estimators ( mi l, and the dashed lines are Monte Carlo 
values using pure pseudo-C^ estimators. For the parameters in the right panel, finite-pixelization artifacts impose a "floor" 
on the value of T/S which can be detected, but the value can be made arbitrarily small by increasing the resolution. For the 
parameters in the left panel, we find that increasing the resolution beyond N S ide = 256 results in no significant improvement. 



ing delensing algorithms. We have crudely studied this 
regime by assuming perfect delensing and zero beam size; 
under these assumptions, we find that the estimators 
alone impose no limit to the value of (T/S) which can 
be detected. 



VIII. CONCLUSION 

In this paper, we have defined pure pseudo-C^ estima- 
tors, which have the property that the estimated B-mode 
power receives no contribution from E-modes, even on a 
cut sky. As a consequence, E-mode signal power does not 
contribute to the variance of the B-mode estimators. The 
pixel weight functions for these estimators have spin-1 
and spin-2 components l|25|) . which can be obtained from 
the spin-0 component by covariant differentiation if its 
form is analytical, or finite differencing (Appendix lB|l in 
general. In contrast to the usual pseudo-CV formalism, 
the spin-0 weight function must have a finite apodiza- 
tion near the boundary, and obey Dirichlct and Neumann 
boundary conditions. We give an algorithm (Appendix 
|P)l for computing the pseudo-CV transfer matrix which is 
used for debiasing, and show that its computational cost 
is C(^maa;)i which is the same as the cost of evaluating 
the estimators. 

We have studied these estimators in detail for mock 
surveys on a 13° spherical cap, using both homogeneous 



noise, and inhomogeneous noise of a specific form 13511 . 
In both cases, we found that for B-mode power spectrum 
estimation, pure pseudo-CV estimators performed slightly 
better than the unmodified versions at noise level 20 /iK- 
arcmin, and much better at noise level 10 /xK-arcmin. 
In the homogeneous case, we considered a wide range 
of noise levels, and showed that the pure estimators are 
~ 80% of optimal, defined by the degradation in (T/ 
which can be detected. The 80% figure is obtained in 
two limiting cases: assuming no separation of the grav- 
ity wave and lensing signals, and perfect separation. In 
the latter case, there is no limit, imposed by the estima- 
tors alone, to the value of [T/S) which can be detected. 
In constrast, using unmodified pseudo-C^ estimators, we 
have found that the smallest signal which can be detected 
in a patch with this geometry is (T/S)i a = 0.042, in 
agreement with Challinor & Chon [8j. 

In future work, we plan to investigate two basic is- 
sues left unaddressed by this paper. First, we will study 
the prospect of generating optimal weight functions for a 
survey algorithmically, starting from Monte Carlo simu- 
lations of the noise, rather than choosing them hcuristi- 
cally. The need for such a procedure is highlighted by the 
inhomogeneous noise model considered in WII in which 
we were forced to use four weight functions with different 
levels of apodization, in order to achieve near-optimal E- 
mode estimation in all £ bands. Second, we will study 
pure pseudo-CV estimators for more realistic noise mod- 
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els, including models which include correlated noise. 
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APPENDIX A: SPIN-1 SPHERICAL HARMONIC 
TRANSFORMS 



As discussed at the end of iUVBI our algorithm for 
fast evaluation of Bg m requires an implementation of 
fast spherical harmonic transforms for vector fields. This 
will also be an ingredient in the algorithm, to be pre- 
sented in Appendix IDl for computing the transfer matri- 
ces K^ e f urc . In this appendix, we write down the recur- 
sion and initial conditions which are needed: 



( C o S e)Yg l -^ 7 ^ T Y l 



(Al) 



. J- 1, 



-(i<W)^_i, m 



C\b 



(Xu )a — £ab{Y u ) 



(-If /(2* + l)! 



2 l t\ y 4ir(£ + 1) 

(sh/- 1 0)e^ [(cosd)x a + iy a ] 
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Here, ( x X lm ) = (l/l) ^/(P - 1)(£ 2 - m*)/{±P - 1). 

An alternate approach to implementing fast spin-1 
transforms appears in [22^ . 



APPENDIX B: FINITE DIFFERENCING IN AN 
IRREGULAR SPHERICAL PIXELIZATION 

Our definition of Bg rn (|24|l requires a prescription for 
computing the spin-1 and spin-2 weights W a , W a b jH^I 
from W(x). If W(x) is of known analytical form, then a 
trivial "prescription" consists of carrying out the covari- 
ant derivatives by hand; otherwise, a finite differencing 
scheme must be used. Since this is nontrivial in an ir- 
regular spherical pixelization such as Healpix, we present 
the details of one such scheme here. 

At each pixel center x, we define a matrix M a b(x) by 



M ab (x) 



(x'-x)^(x'-x)- 



(Bl) 



Throughout this appendix, the notation denotes a 
sum over pixel neighbors x' to the pixel x, and (x' — x) 
denotes projection of the 3-vector (x' — x) into the plane 
perpendicular to the unit vector x. After the projection, 
(x' — ~x) is a tangent vector at x, and M is a rank-2 
symmetric tensor, as the index notation suggests. 

Now, if M^(x) is any pixelized map, we define its "finite 
difference" gradient by 



7 FD 



dcf 



H/(x) = M Q V(x) x 



(B2) 



Yf(x')-/(x))(x'-x 



±b 



For purposes of this paper, we need the following exten- 
sion. If v a (x) is a finite pixelized (tangent) vector field, 
we define its "finite difference" covariant derivative by 



iFD 



def 



v 6 (x) = M- A (x)x 



(B3) 



(Px^xV(x') - U(X)) 6 (X' - X)" 



where 'Px'^xv(x') denotes the tangent vector at x ob- 
tained by parallel translating u(x') along the great cir- 
cle arc connecting x and x'. This parallel translation is 
needed to covariantly difference tangent vectors at dis- 
tinct points x, x'. If one were to simply difference their 
components in the 6, ip coordinate system instead, then 
the resulting finite differencing scheme would suffer from 
coordinate artifacts near the poles. 

Using this finite differencing scheme, one takes the 
spin-1 weig ht W a in Eq. © to be V^ D W, and the 
spin-2 weight W a b to be the traceless symmetric part 
of V^ D {V^ D W). We have tested this prescription for 
the case of cosine apodization, by comparing the pure 
pseudo-C^ estimators which result from finite differenc- 
ing and analytic differentiation of the cosine weight func- 
tion. For a fixed random CMB realization, we find that 
the difference between the two versions is negligible. 



APPENDIX C: CORRELATION FUNCTIONS 
FOR FIELDS OF ARBITRARY SPIN 

For the transfer matrix calculations in Appendix [D] 
we need expressions for correlation functions in terms 
of power spectra, for all combinations of spin-0, spin- 
1, and spin-2 fields. We present these results in a form 
which generalizes to higher spins as well. A spin s field 
(— oo < s < oo) is a function ( s f) whose value at x 
depends on a choice of orthonormal basis vectors {ei, £2} 
at x. Under the right-handed rotation 



= (cos#)ei + (sin#)e2 
e' 2 = — (sin 6)ei + (cos 9)e.2 



(CI) 



( s f) must transform as ( s /)' = e ls6 ( s f)- There is a 
spin-raising operator d and a spin-lowering operator d 
which transform a spin s field into fields of spin (s + 1) 
and (s— 1) respectively. In the frame {ei, 62} = {x a , y a }, 
these operators take the form 



6( s f) = -sin s (0) 



d_ 

09 



icsc(d)^)si,r -(«)(,/) 



sin- s (<9) 



m 



icsc (e)—) S m s (e)( s f) 



An orthonormal basis for spin s fields is given by the spin 
harmonics s Yf m [T(il.l25|]. which are given for s > by 



'(I -a)! 



(C2) 



8*£m 



These are nonzero only for £ > \s\. 

Correlation functions between fields of spin s, s' can 
be expressed in terms of the following function, which 
generalizes the Legendre polynomial Pg (x • x' ) in the case 
s = s' = 0. 

P/ S '(X-X') = f Yl (sYe m (x) { e u e 2} = { X,Y}Y 

m= — t 

x( s /Y^ m (x') { 6 1 ,6 2 }= { x',K'}) (C3) 

We emphasize that the spin harmonics on the right-hand 
side are evaluated in the "two-point" frame {ei,e2} = 
{X,Y}, {ei,e 2 } = {X',Y'}, not the frame {e l7 e 2 } = 
{x a ,Ua}- With this choice of frame, the right-hand side 
is rotationally invariant and therefore depends only on 
the separation (x • x'), as implied on the left-hand side. 

It can be shown [H Eq. (3.17)] that P/ s ' is related 
to the reduced Wigner D-function d l ss , by Pf 8 (cos 6) = 
(—l) s dj nl (Q) Using standard results on Wigner D- 
functions |36(, one obtains the spin-label symmetries, or- 
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thogonality relation, and product rule: 

Pt\z)=Pi' s {z) = (~l) l+s Pr s '(~z) (C4) 



dzPH (z)P/ 2 s (z) 



2L +1 



^2 



(C5) 



p; i is ' i (z)p; 2 24 (z) = x ( C6 ) 



E(^ + i) 



-si -s 2 si + s 2 



s' 2 s[ + s' 2 



p, 



Sl+S2,S,+,S 2 



Likewise, the following recursion and initial conditions 
may be used to evaluate P/ s . Using the spin-label sym- 
metries <|C4|1 . it suffices to give initial conditions for the 
recusion assuming s > \s'\. 



£{£ + !) 

-pY'p!-i{< 



Pt(z) (C7) 



PI" (*) 



(-l) s 



2.s 



1/2 



(l + z) (s+;5 ' )/2 (C8) 



x(l-z)( s - s ')/ 2 (s>\s'\ 



In the recursion l|C7|l . we have defined p\ s — 
^(£ 2 - s 2 )(£ 2 - s' 2 )/£. 

The recursion (|C7jl with accompanying initial condi- 
tion 1C8|I is the main result of this appendix. We have 
used the language of spin-s fields, which permits the re- 
sult to be stated in a uniform way for all spins. Since 
the transfer matrix calculations of Appendix [D] use the 
(equivalent) tensor language in spins < 2, the rest of the 
appendix is devoted to translating between the two. As a 
byproduct, we will obtain expressions for the correlation 
functions between all combinations of fields with spins 
s, s' < 2. In the case where both spins are nonzero, it 
will be convenient to use the linear combinations 



P/ s +(-l) s P, 



s' n>s,-s' 



(C9) 



psi 

TJSS' r l 



-l) s P\ 



instead of P/ s , P/'~ s . 

The basic observation which relates the tensor and 
spin-s formalisms is that the frame-dependent vectors 



m = (ei + ie 2 )/2 



m=(ei-te 2 )/2 (CIO) 



have spins 1 and -1 respectively. (This notation follows 
Okamoto & Hu but our normalization differs by a 
factor 1/V2.) In terms of these, the vector and tensor 
harmonics Y^, can be written 



1 (lm)a 
Y E 



(_iYf m )m a - ( 1 Y lm )m. a (Cll) 
-(-2^m)m b m c - ( 2 Y em )m b m c . 



These are parity-even combinations of {m a ,m a } and 
{m b m c ,m fc m c }. 

The following sums can be evaluated using iRTTTl and 
the definition of Pf s ': 



£y/ m (x)y, ro (xO = ^pf(z) 



4tt 



(C12) 



E Yfyn ( x )^(? m )a' ( x ' ) 



2£+l 

47T 



E Yi m (x)Y ( f m)b , c , (x') — — P® 2 (z)Q' b , c , 



Air 

m 

E Y (?m)a( X ) Y (? m )a'( X ') = ^ X 
m 

(Ql 1 (z)X a X' a ,+Rl l (z)Y a Y^) 

E Yu^ai^YffmWc' ( X = 4^ X 
rn 

(Ql 2 (z)X a Q' b , c ,+R} 2 (z)Y a Ui, c ,) 

E Yum)bc( x )YumWc' ( x = ^ x 
rn 

{Qf{z)Q bc Q' Vd + Rf{z)Q bc U' Vd ) 

To evaluate analagous sums containing YfiL or Y^, note 
that the 90° vector rotation V a — > — e a bV sends Y^, — > 
A -> F, y ^ (-A); and the 45° tensor rotation 
n Qb -> (-e a c n cfc - e fc c n ca )/2 sends -» Y^, Q - U, 
U — > {—Q). For example, the following sum will be used 
in Appendix IdI 



E Y$m)a( X )Y(J> m )bc( X -') — 4- X (C13) 
m 

{Qf{z)X a U' h , d -Rf{z)Y a Q' b , d ) 



It is obtained from the fifth equation in (|C12fl by applying 
a 45° rotation to the indices be, and no rotation to the 
index a. 

As an application, using sums of the form (|C12|I . one 
can write down all correlation functions between Gaus- 
sian fields with spins < 2 which arise from parity-even 
power spectra: 



<^> = E ^<rm 



(C14) 



(xx) = Y j 2 -^{c? g Q\ 1 {z) + c? c rY{z) ] 



(TX) 



E 



4tt 



£ 

(<9Q) - £^(cf B Q2 2 (*) + cf fl ^ 2 (*)) 
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( x Q) = E^(^^ 12 w + ^ B ^ 12 w) 

i 

Here, correlation functions on the left-hand side are de- 
fined in the "two-point" basis from SjH] e.g. (XQ) de- 
notes the correlation between the X a component of a 
spin-1 field and the Q' ab component of a spin-2 field, at 
points x, x' with separation z = (x • x'). 

Correlation functions for parity-odd power spectra are 
given by: 



APPENDIX D: COMPUTING TRANSFER 
MATRICES 

There is still one ingredient missing from our formal- 
ism: an algorithm for efficient calculation of the transfer 
matrices K^ e f ure , which are needed to debias the power 
spectrum estimators (Eq. (|10|l % ). We have found it more 
convenient to compute transfer matrices using correla- 
tion functions than in harmonic space. Since this ap- 
proach has not been used previously in the literature, 
we first illustrate the method for the simpler case of 
temperature pseudo-Q estimators f flP 1|) . before treat- 
ing the case of pure polarization estimators f flP 2(1 . In 
HP 31 we address implementational issues and show that 
the computational cost of computing the transfer matri- 
ces is 0(£ 3 max ). 



(TY) 
(TU) 
(XY) 
(XU) 
(YQ) 
(QU) 



E 



(C15) 



E 

E 

i 

E 

i 

E 



47T 

s 

2^+1 

47T 

21+ 1 

47T 

2£+l. 

47T 

2£+l 

47T 



C^{Qf{z)-R?{z)) 



(cr B Q] 2 (z)-crR^(z)) 



CE u 12r 



{C? E QY{z)-CrR?(z)) 



GB u 12i 



Cf B {Qf(z)~Rf{z)) 



For spins and 2, these results have already appeared in 
the literature |35|,|4(|, with the following explicit formulas 
for P® 2 , Q 22 , Ft 22 used in place of the recursion relation 
dC7lt : 



Qf(z) 

Rf(z) 



= 2 



= -4 



^-i^-(l J ? + ¥)^) 
[(£-l)^+l)(^ + 2)] 1 /2 

(l+2)z p2 , s ( i-4 , t{t-l) \ p 2(„\ 

(e-i)e(e + i)(£ + 2) 

(£ + 2)P 2 _ 1 (z)~(£-l)zP 2 (z) 



(£-!)£(£ +!)(£ + 2) (1 - z 2 ) 



(C16) 



(See also [t| for some spin-1 equations.) The advantage 
of the present treatment is that it generalizes straightfor- 
wardly to all spins, including spins higher than 2, and al- 
lows use of a single uniform recursion relation (|C7|) rather 
than many explicit formulas of type (|C16|I . Our approach 
is similar to Ng & Liu |26( , who compute correlation func- 
tions covariantly for spins and 2. 



1. Transfer matrices for temperature pseudo-CVs 

In the temperature-only case, pseudo-C^ estimators are 
constructed in a way which is completely analagous to the 
polarization case (Eqs. l|8)l-i|lipl. We briefly summarize 
the necessary definitions; details can be found in Hansen 
et al. [3. 

Pseudo multipoles and power spectra are defined by 



T im = ^T(x)W(x)Y^ m (x)* 

X 

1 _ _ 

E/Ti* rri 



a 



2£- 



(Pl) 
(P2) 



m=-e 



Assuming no noise bias, the (£ m ax)-by-(£ max ) transfer 
matrix Ku> is defined by 



TT\ 



E Ku'Cj' 



TT 



(P3) 



The purpose of this subsection is to compute Kw . We 
first establish a key result: evaluation of a sum of the 
form 

S = ^V(x)F(x • x')W(x'), (P4) 

xx' 

where F(z) is a polynomial in z of degree < £ m ax- 
Expanding F in Legendre polynomials, F(z) — 
Y2e<e F(Pi(z), we evaluate the sum as follows: 

S = E W{x)W{itf)F t P e (x-yf) (D5) 



E 

xx'm 



47T 

27TT 



w(x)w(x')F,y, m (x)*y, m (x') 



E w tm wi 

e<e max 



Air 
'21+1 



F 
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It will be useful to rewrite this as an integral involving 
F(z) instead of a sum involving Ff. 



S = 8tt z 



dz( ww (z)F{z), 



(D6) 



We have now arrived at our desired expression for Kggi. 
We defer the issue of efficient evaluation of Eq. (|D8|I to 
9D 31 For now, we note that Eq. l|D8Jl is equivalent to the 
form commonly seen in the literature [e.g. ll5l Eq. A31], 



where 



-WW 



(D7) 



can be thought of as an estimator of the correlation func- 
tion of W(x), pretending that W(x) is a Gaussian field 
on the full sky which is zero (by coincidence!) outside 
the survey region. 

Armed with the key result 1|D6|1 , the transfer matrix is 
calculated as follows. By definition, the matrix element 
Ku> is the expectation value of Cj T (Eq. JEJ)), given the 
signal covariance (T(x)T(x')) = (2f + 1)/(4tt)P^ (x ■ x'). 
Therefore, 



K, 



2f 



2tt 



1— 53 M/(x)^(x')r, m (x)y,„ l (x') , 



2C 



4tt 



(D8) 



1 dz C WW (z)Pi(z) ri '' 



47T 



Ku., = 



5J |WW| 2 

t"m" \ 




(D9) 



by virtue of the definition 1|D7|) of <^ ww (z) and the iden- 
tity 



dzPi{z)P v {z)P l »{z) = 2 



i e i" 





(D10) 



2. The transfer matrix KfX ure 



Computing the transfer matrix K^f ure defined in Eq. I|29|l is completely analagous to the temperature-only case, 
but the bookkeeping is considerably more complicated owing to the presence of spin-0, spin-1, and spin-2 weights 
W(x), W a (x), Wbc(x). We will need to use the results in Appendix ICl in which it is shown how to compute correlation 
functions between fields of these spins from all parity-even power spectra. 

Throughout this subsection, we will use the abbreviated notations 

Y lm d ^ Ytm (x) Y' lm ^ Y lm (x') z d ^ (x • x') 

in equations which contain a pair of points x, x'. 

We first consider an analog of the "key result" Eq. 1D4|I , in which F is replaced by a 3-by-3 matrix which contracts 
all combinations of spin-0, spin-1, and spin-2 weights: 



where 



F(x,x') 



S = J2(w(x) W a {x) W bc {x) ) F(x,x') 



/ W(x) 

VF a '(x') 
^ W b ' c '(x') 



(Dll) 



f WW (z) f wx (z)x' a , V WQ {z)Q'v c . 

f wx {z)X a f xx (z)X a X' a , + f Y {z)Y a Y' a , 2f x Q(z)X a Q' b , c , + 2f YU {z)Y a U' b , c , 
2f w( 3(z)Q bc 2f x Q(z)Q bc X' a ,+2f YU (z)U bc Y>, 4fQQ(z)Q bc Q' b , c , + Af uu (z)U bc Ui, c , 



Note the bi-tensor structure of F: each matrix entry carries a set of tensor indices which contract with the weight 
functions in (|D 1 1|> so that S is a scalar. 

We make the assumption that F is built from parity-even power spectra which are zero for multipoles larger than 
P-max- (This is the analog of the assumption, from the preceding subsection, that F(z) is a polynomial of degree 
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< £ max .) To be precise, let Ff G , Ff E \ Ff E , Ff E , Ff c , Ff B , Ff B be arbitrary power spectra. Note that 

we denote the spin-0 "field" by W, and consider only parity-even combinations of fields. Then F is assumed to be of 
the form, 



F(x,x') 



E 



/ tpWWv* yi 
I r l I lm 1 lra 



Kir, 



rpWGy* viG 
r l I im I (tm)a' 



rpWWy* yiE 
r i 1 im I (lm)b'c' 



pWGyG* vi pGGvG* v/G pGEvG* v/E 

r l I (im)a I im r t 1 (f.m)a 1 (lm)a' r l 1 (tm)a L {lm)V d 



rlE 



i jpWEvE* vi rpGEvE* viG rpEEvE* yi 

\ r i 1 (im)bc 1 im r l 1 (im^c 1 \lm)a' r l 1 (£m)fac 2 (£m)b' c< 



+ 



( 

n pCCvC* y'C 

U r i 1 (Im)a 1 (tm)a' 

\ pCBvB* y'C 

\ U e l 1 (im^c 1 (im)a' 





pCByC* yiB 
r l 1 {lm)a (tm)V d 
TpBBvB* yiB 
C l 1 {lm)bc {lm)b' d 



(D12) 



The significance of the sum (|D11|I . with F-matrix of the form (|D12I) . will appear shortly. For now, we forge ahead 
and evaluate S, using the method of the preceding subsection. First, plugging the form (|D12Jl for F into the definition 
(IDTT) of S, we get: 



E ( W ^ 



wL wf m 



VV Im VV im 



/ pWW pWG pWE 
pGG pGE 



WG 
WE 



F? E F, 



EE 



( WT 

im 

w9* 

im 

\ W E * 



pCC pCB 
pBB pBB 



(D13) 



This expression for S is the analog of Eq. (|D5(I from the preceding subsection. Our goal is to get an expression 
analagous to Eq. (HSU, in which the functions {f ww (z), f wx (z), ...}. appear directly instead of the power spectra 
{Fj vw , F(* i ■ ■ ■}■ We first define correlation functions which are constructed from the weight functions in the same 
way that one constructs correlation functions from parity-even power spectra (Eq. (|C14I) ): 

(D14) 







(WW 


im 




im 


+ W^Wg.Rf) 


im 


C YY = 


^^(W^W^Rj 1 
Im 


+ WglW^Qf) 


c WQ = -^Y,( w ^ w fmrn 

im 


(■QQ _ 


Im 


+ W^W^Rf) 


im 


qUU _ 


^Y.iw E *w E m RT 

im 


+ w^wtQf) 


im 



Eq. i|D13|) can then be rewritten as an integral containing the functions {f ww 



S = 8tH 



dz 



(z),f wx (z),...}: 

ww (z)f ww (z) + 2Q WX (z)f wx {z) + 2C WQ (z)f w ^(z) + C XX (z)f xx (z) 
+ C YY (z)f YY (z) + C YY (z)f YY (z) + 2t x <J(z)f x <J(z) 
+ 2( YU (z)f YU (z) + C QQ {z)f QQ {z) + C UU (z)f uu (z) 



(D15) 



To show this, it is easiest to work backwards, substituting Eqs. I|D14(I into Eq. (|D15|I . and ending up with Eq. l|D13j) . 
One writes Q*/ , R s / in terms of P/ s ' (Eq. |H9])). and uses the orthogonality relation for P" (Eq. lC5j) V 

We have now arrived at our desired "key result": an expression for S in terms of the functions {f ww (z), f wx (z), 
. . .}. We now proceed to calculate the transfer matrix element K^ e f wre . (The case of K~[ e f ure will be treated shortly.) 
This matrix element is the expectation value {Cf B,pure ) given the signal covariance 



n de (x)n d ' e '(x')) = (^^) (Rf(z)Q de Q ,d ' e ' + Qf 



! (z)U de U' d ' e '). 



(D16) 
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A short calculation, using only the definition of (j BB -P ure (Eqs. (|24|l . (J2HJ); shows that 

( 



W a (x') 
y W b ' c '(x') 



(D17) 



where the F-matrix is given by 

[Rf{z)Q^' + Qf(z)U de U'^'] x 



F«/(x,x') = 



(lm)de 



N' t T de jY^ ){ ( Y{* m)d , e , N[T dlel jY[f m)r NiT d , & , h , c ,Y' tm ) 



N e T debc Y e * m 

We will evaluate the right-hand side of Eq. (|D17|) using the key result (|D15(1 , but first we simplify the F-matrix using 
identities from Appendix IU1 We show the details for the (1,2) matrix entry, 



F^=A^±1 



iG 

{im)f> 



, [Rf^QteQ'^' +Q¥(z)U de U' d ' e '} Y { ^ )de T d , e , a ,f'Y { 

n{M + L > m=-t 

To do the sum over m, we use identity (|C13|I . obtaining 

F^' 2) = N'^^[Rf{z)Q de Q' d ' e ' +Qf(z)U de U' d ' e ']T d , e , a ,f'(Q 1 e 2 (z)U de X> f , - R?{z)Q de Y' fl ) 
It is then straightforward to do the index contractions and obtain: 



F (i,2) 



16tt 2 v 1 l ' H H ' 

This procedure can be used to simplify the remaining entries of F«/; when the dust has settled, one finds 

21' + 1 



p(2,2) 
f ppi 



P (2,3) 
1 OPI 



p(3,3) 
f ppi 



N, 



Wt2 {QfQf + RfRf) 

£ "l' ^b'c' 



(D18) 

(D19) 
(D20) 
(D21) 



= N, 



2f + 1 
16tt 2 



(Q] <3l' + RYRf)X a X' a , + (R} Qjf + QfRf)Y a Y' a , 



-N e N' e ^±^P^(QfX a Q' b , c! + RfY a U' Vc ,) 
Nf^-P?°(Q?Q bc Q' b , c , + RfU hc U' Vc ) 



We would like to apply the key result (|D15(I with F-matrix given by this form, but first there is an annoying technicality: 
our derivation of Eq. I|D15|I assumed that the F-matrix was built from parity-even power spectra, in the sense that 
Eq. I|D12|) is satisfied. We claim that this is so for the F-matrix in (|D21(I . with £ max =£ + £'. To make this 
statement intuitively plausible, note that the F-matrix in (|D21fl is a rotationally invariant, parity-even object which 
is constructed by multiplying objects of spins £ and I 1 . A formal proof can be given by writing Q\ s , Rf s in terms of 
Pi 3 ' (Eq. jnni), and using the product rule for P/ s ' (Eq. Ipijl l. This can be done one matrix entry at a time; again 
we supply the details only for the (1,2) matrix entry: 

,2f + 1 



r (l,2) 
W 



(x,x') 



N, 



(4tt) 2 
,21' + 1 
2(4tt) 2 



(Q?Qf+R?Rf)X' a , 

{prpF+p}- 2 Pi- 2 )x> a , 



(D22) 



E v 

-e>\<e"<£+e 



l (2£' + l)(2£" + 1) / I (' £"\(l £' £" 
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We have shown that F^; 2 -* is of the form which appears on the RHS of Eq. I|D12|) . with £ max =1 + 1'. 

With this final technicality out of the way, we can use the key result l|D15Jl . with F- matrix given by l|D2T]) . to 
evaluate the right-hand side of Eq. I|D17() . This gives the transfer matrix in the form: 

K+f ure = 2n dz {^-^j (Mz)Qf(z) + B e (z)R 22 (z)) (D23) 

where the functions Ai(z), Bi{z) are defined by 

An = C WW Qf - 2N^ WX Q 12 + 2N e ( w Qp° 2 + N[ 2 (C XX Q} 1 + Q YY R} 1 ) - 2N e N^ XQ P^ + N 2 ( QQ P^ 



B e = ( WW R 22 - 2N' l Q wx R] 2 + n'^^R] 1 + ( YY QV) ~ 2N £ N'( YU P° 1 + N 2 <; UU P° 

To compute K^ ure , one modifies the signal covariance (|D16|I by exchanging Q 22 <-> R 22 ; this modification carries 
through to the end of the calculation and shows: 

K~r re = ^ \\ ^ (j^ 1 ) (Mz)Rf(z) + B e (z)Q 22 (z)) (D24) 

Eqs. I|D23|1 and (|D24|I are the main results of this appendix, and show how to compute the transfer matrices K^ e f ure 
for any choice of weight functions. 

We emphasize that even though the derivation is lengthy, the final result is simple. Each term in the integrals 
(|D23J| . (|D24(I is a product of three functions: a generalized Legendre polynomial of degree £' , a generalized Legendre 
polynomial of degree £, and a correlation function {C (z), Q {z), . . . } which depends only on the pixel weight 
functions l|D14|l . The final form of the transfer matrices is simpler in position space than in harmonic space, where 
the three-way multiplication would be replaced by a sum involving 3j symbols, similar to Eq. (|D9(1 . but with many 
terms. 

Finally, we discuss the case in which multipoles are binned into bandpowers. In this case, the bandpower transfer 
matrix K^ b , is related to K^, by: K^, — P b iK^,Pgi b i , where the matrices P and P define the binning (SjfflJ. For each 
bandpower, we define "binned" versions of the functions Ai{z), B?(z), Q 22 (z), R 22 (z): 



A b {z)=Y J PbiAi{z) 

i 



3. Efficient calculation of the transfer matrix 
integrals 

We have now shown that the transfer matrices K^, can 
be represented in integral form IJD26I) . In this subsection, 
we give an algorithm for evaluating the integrals to ma- 
chine precision, whose running time is 0(^ aa ,). Here, 
£ ma x is the largest value of (£ + £') for which a trans- 
fer matrix K^, must be computed. This is the same 
cost, within a constant factor, of evaluating the estima- 
tors once (SjlVBj). 

The first observation is that, when computing correla- 
tion functions {( ww (z), Q WX (z), . . . } using Eq. (jDWj, 



B b {z)=Y J PbiB t {z) (D25) 

t 

iM*) = £(^W^ 2 (z) 



(D26) 



it is only necessary to sum over multipoles £ < £ ma x- This 
is because, as argued at the end of the preceding section, 
the F-matrix can be written in the form 1D12|) . and by 
Eq. (|D13(1 . only multipoles of the weight functions with 
£ < £max contribute to S. 

The second observation is that the integrands in l|D26|) 
are polynomials of degree < 2£ max . This can be seen 
term- by-term after plugging in the definitions of Ag(z), 
B e (z). For some terms, such as ( ww (z)Qf (z)Qf (z), 
all three factors are polynomials; for others, such as 
( wx (z)Q 12 {z)Q 22 {z), the last is a polynomial and the 
first two are polynomials times \J\ — z 2 . There are no 
terms with an odd number of \f\-z 2 factors. 



In terms of these, the bandpower transfer matrices are given by: 

K+r re = 2tt J 1 ^dz{A b {z)Q b ,{z) + B b {z)R b ,{z)) 
K- b f ure = 2vr J^dz(A b (z)R b ,(z) + B b (z)Q b ,(z)) 
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Because of this second observation, the integrals can be 
done exactly using Gauss-Legendre quadrature §4.5] 
with (£ max + 1) points. Our algorithm for evaluating the 
integrals is therefore given as follows. First, we com- 
pute spherical harmonic transforms {W^ TO , W^, W^ ml 
^imi °f the weight functions. Second, we com- 
pute the correlation functions {( ww (z), ( wx (z), ...} 
at each of the (£ m ax + 1) quadrature points, using Eq. 
(|D14|I . Third, we compute the AN^and values {At(z), 
Bb(z), Qb( z ), Rb( z )} at each quadrature point z, us- 
ing (|D25(I . Fourth, we loop over bands b, &', comput- 
ing K b il b f ure by doing the integrals in (|D26fl by Gauss- 
Legendre quadrature. 

Let us consider the running time of each of these 
stages. The first stage is 0(1^^), using fast spherical 
harmonic transforms. The second stage can be done in 
time 0(£ 2 nax ) by first summing over m in Eq. (E2J) 
at fixed £ and then summing over £ at fixed z. The 
third stage is 0{Nf ]an d£^ nax ) 1 evaluating the generalized 
Legendre polynomials by recursion (|G7|I . The fourth 
stage is 0{N 2 and £ rnax ) since Gauss-Legendre quadrature 
is 0(£ max ) once everything has been precomputed at the 
quadrature points. Putting this together, and noting 
that Nband < £max, the computational cost of computing 
the transfer matrices K^ e f ure is 0(£^ nax ). 



APPENDIX E: COMPUTING NOISE BIAS 

In addition to the transfer matrix K^ e f ure , the estima- 
tors presented in this paper also require computing noise 



bias terms N EE , N BB ' pure (Eq. JHJ}). For real exper- 
iments, which include such complications as 1// noise, 
noise bias for all types of pseudo-CV estimators must be 
computed by Monte Carlo. Indeed, a practical advantage 
of the pseudo-CV framework is that unbiased estimators 
can be constructed given only Monte Carlo simulations 
of the noise; no other representation of the noise covari- 
ance is required. However, for theoretical studies, it is 
convenient to have an exact formula for the noise bias in 
simple cases. In this appendix, we consider noise which is 
uncorrelated between pixels, and isotropic in each pixel, 
but not necessarily homogeneous: 



(Q(x)Q(x)) = <C/(x)[/(x)) = a(x) 2 <5 x 



(El) 



Here, we represent the noise by its per-pixel RMS tem- 
perature cr(x) (i.e., units [iK rather than /iK-arcmin). 

With noise covariance given by fElJ, the noise bias 
defined in l(2l)|) is given by 



N, 



N EE 



BB,pure 



E 

X 

E 4 -j 



^ 2 (x) 

4-7T 

a 2 (x) 



TT(x) 2 (E2) 
W(x) 2 + N' e 2 W a (x)W a (x.) 



+2N 2 W bc {*)W bc {*) 

This is derived starting from the definitions of Ei m , Bi 
(Eqs. ©, |23I)) using identity (|C12j) from Appendix 



APPENDIX F: FISHER MATRIX EVALUATION WITH AZIMUTHAL SYMMETRY 

In this appendix, we present the details of our method for fast exact evaluation of the Fisher matrix 1)31(1 , in the 
case of inhomogcneous, but azimuthally symmetric, noise. A similar method, in the context of CMB temperature, 
appeared in |2jfl . Since we only consider uncorrelated noise in this paper, the noise covariance can be written 

(Q(x), Q(x')> = (C/(x), C/(x')) = v(0) 2 S i2) (x - *'), (Fl) 

where "q(9) is arbitrary. If we change variables from {Q, U} to = (Q±iU), and Fourier transform in the azimuthal 
coordinate tp, by defining 

r 2ir 



H±(#)= / dc P {Q±tU}(6,c P y m * (F2) 
Jo 

then the noise covariance (|F1|) is still diagonal (in both m and 9): 



<n+*(0)II+ {9')) = (U-*(9)U-,(9')) = 4*^5(9 - 0')<W (F3) 

The point of this change of variables is that the signal covariance is also diagonal in m (but still dense in 9) . Using 
the results of j^l), one can show that the signal covariance is given by: 



n-*(#) 



(n+,(0O u^(9')))=A^s mm , x (F4) 



V / (Cf E + Cf B )UY lm {e^)){-2Y lm {9'M) [Cf E -Cf B )UY em (8,0))(2YU0',0)) 
4" V (Cf E - Cf B )( 2 Y em (9,0))UY hn (9',0)) (C EE + C BB )( 2 Y em (9,0))( 2 Y tm (0',0)) 
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Eqs. (|F3|) and l|F4(l can be used to efficiently calculate the Fisher matrix 

F w = iTr(S b (S + N)"% (S + N)" 1 ) (F5) 

For concrete calculation, the continuous coordinate is replaced with a set of N equally spaced values 0\ , 2 , ■ ■ ■ , 
On, with spacing AO. The (2JV)-by-(2iV) noise and signal covariances of the variables {11+ (Oi), • • • , n+(0jv), II~ (#1), 
• • • , Il~ (0jv)} are then given by 

N (m) _ il / 3§$J*« \ (F6) 

^sj (F6) 

g(m ) = 4?r 2 ( {Cf E + C BB ) - 2 Y em (0 h 0) - 2 Y lm {0 h 0) (Cf s - Cf B ) _ 2 ^ m (^, 0) 2 ^ m (0j, 0) \ 
77 Y \ {Cf E -Cf B ) 2 Y im {0 i ,Qi)^Y im {0 j ,Q) (Cf E + C BB ) 2 Y em (0. l ,0) 2 Y em (0 J ,0) J 

The Fisher matrix (|F5|) is given by summing over m and tracing over 6: 

*W = ^ £ TrCS^CS^+^^S^CS^+N)- 1 ) (F7) 

rri--m ma x 

h rrimax 

= iTr(S[ 0) (S( 0) + N)- 1 Si? ) (S^ 0) + N)- 1 ) + £ Tr(s£ m) (S<, m) + N)- 1 S^ ) (S^ ro) + N)" 1 ) 

m— 1 

where "Tr" is a (2iV)-by-(2iV) trace. Evaluating the spin harmonics by recursion in £, the computational cost of 
each Fisher matrix element is 0(Ngm max ), versus 0(N^ ix ) for brute force calculation without exploiting azimuthal 
symmetry. (For the mock surveys studied in this paper, with Healpix resolution N S id e — 256, one has Ng ~ 70 and 
N pix ~ 10 4 .) 



